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Abstract 

We study hydrogen-bond dynamics in liquid water at low temperatures using 
molecular dynamics simulations. We find that bond lifetime ("fast dynam- 
ics") has Arrhenius temperature dependence. We also calculate the bond 
correlation function and find that the correlation time ("slow dynamics") 
shows power-law behavior. This power-law behavior, as well as the decay of 
the bond correlations, is consistent with the predictions of the mode-coupling 
theory. The correlation time at the lowest temperature studied shows devia- 
tion from power-law behavior that suggests continuity of dynamic functions 
between the liquid and glassy states of water at low pressure. 
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Both experiments and simulations |5|-|7j of water have focused on understanding 

various aspects of hydrogen bond dynamics, such as the average bond lifetime t hb and the 
structural relaxation time tr. Experiments show that tub follows a simple Arrhenius law, 
but that characteristic relaxation times commonly obey power laws on supercooling 
Simulations are particularly useful for investigating water in the supercooled regime since 
nucleation does not occur on the time scale of the simulations and access to quantitative 
hydrogen bond information is available. Accordingly, we carry out simulations of normal 
and supercooled water and find Arrhenius behavior of thb, as expected from experimental 
results fHH; the activation energy associated with tub provides a simple test of the criteria 
used in the simulations to define a bond. Further, we find power-law behavior of tr which 
breaks down at the lowest temperature. This non-trivial growth of tr is consistent with a 
continuous transition of the supercooled liquid to the glassy state. 

To elucidate how simple temperature dependence of tub [Fig. [| can lead to complex 
dependence of tr [Fig. we perform lengthy molecular-dynamics simulations (up to 70 ns, 
over 7 orders of magnitude) at seven temperatures between 200 K and 350 K ||. We 
use the extended simple point charge (SPC/E) potential |TIJ for water and calculate the 
bond dynamics by considering two definitions of an intact hydrogen bond: (i) an energetic 
definition ||, which considers two molecules to be bonded if their oxygen-oxygen separation 
is less than 3.5 A and their interaction energy is less than a threshold energy Ehb, an d (ii) 
a geometric definition , which uses the same distance criterion but no energetic condition, 
instead requiring that the 0-H...0 angle between two molecules must be less than a threshold 
angle 9 HB . 

We calculate tub using both bond definitions over the entire temperature range simu- 
lated and find Arrhenius behavior of tub [Fig- fl- Measurements of tub using depolarized 



light scattering techniques |2|,[| find Arrhenius behavior [12|]. Previous simulations did not 
find Arrhenius behavior of tub @- The activation energy Ea associated with tub has been 
interpreted as the energy required to break a bond via librational motion, a "fast" mo- 
tion . Comparison of experimental and simulated values of Ea provides a primitive test 



of the bonding criteria in our simulations; we obtain reasonable agreement between experi- 
mental and simulated values of Ea using thresholds of Ehb = — 10 kJ/mol for the energetic 
definition and 6hb = 30° 0] for the geometric definition. We find better quantitative agree- 
ment with experiments for tub values obtained from the geometric definition than for tub 
values obtained from energetic definition - possibly because the geometric bond definition, 
like the depolarized light scattering experiments, is highly sensitive to the linearity of the 
bond. We also calculate tub using the thresholds Ehb = kJ/mol f|]f| and Qhb = 35° and 
find Arrhenius behavior, but with Ea roughly 30% smaller for the energetic definition, and 
roughly 10% smaller for the geometric definition. 

The quantity tub is one characteristic time - the mean - of the distribution of bond 
lifetimes P{t), which measures the probability that an initially bonded pair remains bonded 
at all times up to time t, and breaks at t. P(t) is obtained from simulations by building 
a histogram of the bond lifetimes for each configuration, and may be related to tub by 
Tub — Io° tP(t)dt ||16|| . For both bond definitions and all seven temperatures simulated, we 
calculate the behavior of P(t) [Fig. 0(a)] . We observe neither power-law nor exponential 
behavior for either bond definition, unlike previous calculations using the ST2 potential 
which revealed power-law dependence of P(t) |]]. The difference in P(t) is not surprising, 
since the previous study considered a different threshold value and also followed a path near 
a liquid-liquid critical point [|l7],[il|. The behavior of P(t) for the two bond definitions is 
different, likely caused by differences in sensitivity of the definitions to librational motion. 

We study tr ("slow" dynamics) for an initially bonded pair by calculating the bond 
correlation function c(t), the probability that a randomly chosen pair of molecules is bonded 
at time t provided that the pair was bonded at t = (independent of possible breaking in 
the interim time). To calculate c(t), we define 

c(t) = (h(0)h(t))/(h 2 ), (1) 

where h(t) is a binary function for each pair of molecules {i, j}, and h(t) = 1 if molecules 
{i, j} are bonded at time t and h(t) = if {i, j} are not bonded at time t 0. The angular 



brackets denote an average of all pairs {i, j} and starting times. We define tr by c{tr) = 
e _1 . Short time fluctuations and choice of bond definition do not strongly affect the long- 
time behavior of the "history-independent" quantity c(t) (which does not depend on the 
continuous presence of a bond), but such fluctuations cause both qualitative and quantitative 
differences in the long-time behavior of the "history-dependent" P(t) (which does depend 
on the continuous presence of a bond). 

We can interpret the behavior of tr in terms of the mode-coupling theory (MCT) for a 
supercooled liquid approaching a glass transition Hlj|l9|. In accordance with MCT, we find 
- independent of the two bond definitions considered - power- law growth for T > 210 K 

t r ~(T-T c )-\ (2a) 

with 7 = 2.7 ± 0.1 and T c = 197.5 ± 1.0 K, approximately 50 K less than the temperature 
of maximum density Tmd of SPC/E [Fig. 0(a)] p0| |. Fits of experimental relaxation times 



to Eq. ( |2a|) also find T c at a temperature 50 K less than the Tmd of water ||21|| . In MCT, 
T c is the temperature of the ideal kinetic glass transition and is larger than the glass tran- 
sition temperature T g determined from, e.g., viscosity or relaxation time measurements. At 
T = 200 K, tr is smaller than would be estimated by Eq. (|2a|), most likely because MCT 
does not account for activated processes which aid diffusion and reduce relaxation times at 



low temperatures f22 |. Typically, these activated processes become important near T c , as 
observed here. 

Our simulation results for tr can also be fit by the Vogel-Fulcher-Tammann (VFT) 
form for T > 210 K 

T R ~e A /( T - T °\ (2b) 

with T = 160 K [Fig. 0(b)]. In the entropy theory of the glass transition P|J2^|, T is 
associated with the Kauzmann temperature, the temperature where the extrapolated en- 
tropy of the supercooled liquid approaches the entropy of the solid. For typical liquids we 
expect T < T g , so estimating T and T c provides lower and upper bounds for T g . For water, 



however, fits to Eq. ( pliiD of experimental data (which are far above T g ) yield T > T g [|TJ; 
hence we do not consider our T value to be a lower bound for T g of the SPC/E model. 
T g is defined experimentally as the temperature where the viscosity reaches 10 12 Pa-s or 
t = 100 s. Experiments near T g often show a crossover from VFT behavior to "normal" 
Arrhenius behavior [|IJ. While our simulations are still relatively far from T g (based on the 
value of tr), a naive extrapolation assuming that temperatures T = 210 K and 200 K follow 
Arrhenius behavior yields T g « 105 K ||24j| . Applying the same shift to T g as is observed 



for the Tmd of SPC/E relative to water (specifically 35 K), the speculated T g of SPC/E is 
consistent with experimental measurements of T g « 140 K in water |]J . 

A continuous crossover to Arrhenius behavior in water might account for the fact that 
T g < T when experimental data are fit to Eq. (|2b[) . VFT or power law behavior ("fragile") 
for T > 220 K changing to Arrhenius behavior ("strong") for T < 220 K could smoothly 
connect the structural relaxation times in the liquid with those of the glass. A fragile-to- 



strong transition in water has been previously suggested ||26|| , and recent experimental results 
for the diffusion constant (which scales as r^ 1 ) at temperatures closer to T g may help to 
determine if such a transition occurs lEH . 



The reactive flux, defined by the derivative 

k(t) = -dc(t)/dt, (3a) 

measures the effective decay rate of an initial set of hydrogen bonds. At temperature T = 
300 K, non-exponential decay of k(t) was found for the closely-related SPC model using 
the geometric bond definition 0. Our calculations of k(t) [Fig. 0(b)] reveal a power-law 
region k(t) ~ for T > 250 K for both bond definitions, with an exponent £ = 0.5 ±0.1. 
The range of the power-law region increases from about one decade at 350 K to about two 
decades at 250 K. 

The value of ( can be interpreted using MCT, which predicts that c(t) decays from a 
plateau value c p with power-law dependence 

c p - c(t) ~ t b (3b) 



in the range where k(t) appears to be power-law [p],|19[. From Eq. (|3b"D , k(t) ~ t , so 
b = 1 — We find b = 0.5 ± 0.1, consistent with previous work ||11|| , further suggesting that 
the bond behavior is consistent with MCT predictions for a glass transition. For T < 250 K, 
the decay of k(t) appears to be neither power-law nor exponential. 

We note that even though the functional form of thb does not appear to be strongly 
dependent on bond definition, P(t) is different for the two definitions — suggesting that 
P(t) may not be the best function for studying bond dynamics. In contrast k(t) — which 
includes bond reformation — appears to be largely independent of the bond definition at long 



times [28|| . 

The Arrhenius behavior and short time scale of thb indicates that librational motion 
is a thermally activated process that is largely independent of the structural slowing. The 
non-Arrhenius behavior of tr and the non-exponential relaxation of k(t) are consistent with 
the presence of the proposed ideal kinetic glass transition for SPC/E approximately 50 K 



below the Tmd of SPC/E water ]TT[, which coincides with the temperature at which many 
experimental relaxation times appear to diverge [TJI|] . It was hypothesized that the apparent 
singular temperature of liquid water observed experimentally may be identified with the T c 
of MCT jn],^7J]; our results support this hypothesis. They also complement scenarios that 
account for the anomalous thermodynamic behavior, such as the liquid-liquid transition 



hypothesis ||i7| , f2T||30l or the the singularity-free hypothesis 
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FIG. 1. Average bond lifetime thb for the energetic (o) and geometric (x) bond definitions; 
shown for comparison are experimental data (o) for depolarized light scattering Q. We observe 
that all of the data can be fit by Arrhenius behavior tjjb = T oGxp(EA/kT) with approximate 
activation energies: (i) Ea = 8.8 ± 1.0 kJ/mol (energetic definition); (ii) Ea = 9.3 ± 1.0 kJ/mol 
(geometric definition); (hi) Ea = 10.8 ±1.0 kJ/mol (experimental). We scale the temperature of 
the simulation results by T^^ E and temperature of the experimental data by T^ 2 ^ to facilitate 
comparison with with results [13-15]. 
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FIG. 2. Relaxation time tr of the hydrogen bond correlation function c(t) for the energetic (o) 
and the geometric (x) bond definitions, (a) Fit to the scaling form predicted by mode-coupling 
theory (solid line) with T c = 197.5 K. In MCT, T c is the temperature of structural arrest, (b) Fit 
to the VFT form (red line) with To = 160 K. The deviation from both fitting forms we consider 
may indicate a smooth transition relaxations time in supercooled water with that of glassy water. 
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FIG. 3. Time-dependent relaxation of the hydrogen bonds, (a) The bond lifetime distribution 
Pit). The curves can be identified as follows (reading from top to bottom): 200 K, 210 K, 225 K, 
250 K, 275 K, 300 K, and 350 K. Each curve is offset by one decade for clarity, (b) The reactive 
flux k(t). Note that k(t) does not depend on the unbroken presence of a bond, so k(t) decays 
less rapidly than Pit). After a transient period of rapid librational motion up to t ~ 0.3 ps, we 
observe a region of power-law decay for T is above 250 K. For t > 0.3 ps, k(t) is nearly identical for 
both bond definitions, a result which likely arises because both definitions use the same distance 
criterion. We calculate kit) from the numerical derivative of c(t), which is well-approximated by 
a stretched exponential for t > 0.3 ps. We average over all pairs in the system and many initial 
starting times for trajectories ranging in length from 200 ps at 350 K, to 40 ns at 200 K. Note that 
our results for the geometric definition at T = 300 K are consistent with recent calculations for 
the SPC potential (see inset of Fig. 1 of ref. @). 
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